Conditioning 3D Stochastic Channels to Pressure Data

This work discusses the development and implementation of a procedure to condition a stochastic channel to welltest pressure data and well observations of the channel thickness and the depth of the top of the channel. The stochastic channel is defined by a set of geometric random variables (referred to as geometric model parameters) that describe the location, size and shape. Channel and nonchannel permeability and porosity are treated as random variables. These four random variables plus the geometric parameters comprise the complete set of model parameters. Multiple conditional realizations of the geometric parameters and rock properties are generated to evaluate the uncertainty in model parameters and the reduction in uncertainty obtained by conditioning to well-test pressure data. Introduction Georgsen and Omre have presented a stochastic model for river beds (referred to here as channels) within a channel belt and used this model to simulate a fluvial reservoir consisting of multiple channel belts. Except for the use of aspect ratio in place of thickness as one of the Gaussian random fields used to define the size of a stochastic channel, the set of geometric parameters that describe our channels is nearly identical to theirs. Georgsen and Omre also implemented a simulation algorithm to generate conditional realizations of multiple channels within a channel belt where the conditioning data consist of facies observations at each well. Georgsen extended the model of Georgsen and Omre to include four facies, namely a background facies, sheet splay, channel sand and barriers. In the simulation process the last three facies are embedded as geometric, geological objects within the background facies. All facies exist within a rectangular threedimensional box which represents the complete reservoir system. In our work, we simulate a single channel within a rectangular three-dimensional box. However, unlike the two references cited above, the permeability and porosity within and outside the channel are considered to be random variables so the model parameters include both the geometric parameters for the stochastic channel and the rock property fields. Moreover, we generate realizations of the model conditioned to single-phase flow pressure data as well as well observations of the channel. The problem of conditioning a channel to pressure data has been considered previously by Landa and Horne (also see Landa) and Rahon et al. Landa and Horne, however, only considered two-dimensional channels and used a model in which the channel width is constant and channel boundaries are described by trigonometric functions so that the number of channel geometric parameters was very small. Moreover, channel and nonchannel permeability and porosity were assumed to be known. They showed that the sensitivity of pressure to channel parameters can be calculated efficiently using the direct method. (The direct method is sometimes referred to as the sensitivity coefficient method (Yeh) in the hydrology literature and is commonly referred to as the gradient simulator method of Anterion et al. in the petroleum engineering literature.) In the work of Landa and Landa and Horne, channel parameters were estimated by applying the GaussNewton method to minimize an objective function consisting of a sum of squares of data mismatch terms. Various modifications, including Levenberg-Marquardt and modified Cholesky decomposition, were used as necessary to provide regularization. In the single synthetic example presented by Landa, the channel parameters were esti2 CONDITIONING 3D STOCHASTIC CHANNELS TO PRESSURE DATA SPE 56682 mated using two-phase (oil-water) flow production data at multiple wells, drill stem test pressure data and a map of changes in saturation. (Landa assumed the saturation map could be obtained reliably from 4D seismic data.) For this example of Landa, the algorithm required 41 iterations to converge. He attributed the slow convergence to the existence of a local minimum and suggested reweighting of individual data mismatch terms in the objective function to overcome the problem. However, in the probabilistic approach applied here, the objective function to be minimized appears as the exponent of a probability density function and can not be adjusted arbitrarily. For the problems we have considered, we found that convergence properties of the algorithm can sometimes be significantly improved by conditioning to an observation of the channel at the well. We present a procedure to incorporate such information in an efficient way that improves the convergence properties of the algorithm and is consistent with the randomized maximum likelihood method used to generate realizations. Rahon et al. modeled a 3D channel by triangularization of its bounding surfaces. In the three-dimensional example they considered, the channel is actually defined by a triangularized surface with 90 nodes. They consider only the case where the channel lies along a specified line, i.e., the centerline of the channel is a known line. They assume porosity and permeability within and outside the channel are known. Sensitivities are computed by a continuous form of the adjoint method and they estimate channel parameters (the positions of nodes of triangles) by minimizing an objective function that includes the squares of both pressure and pressure derivative data mismatches. In conditioning a geostatistical model to dynamic data, computational efficiency is a key issue. Any method for generating conditional realizations of the model or the maximum a posteriori estimate requires the minimization of an appropriate objective function. Although many references have applied a Gauss-Newton method for optimization, one sometimes finds that it is necessary to damp the change in model parameters at early iterations to avoid slow convergence or convergence to a local minimum which gives an unacceptable match of observed data. This problem has been discussed recently by Wu et al. and is similar to convergence difficulties observed in very early work on automatic history matching; see, for example, Jacquard and Jain and Jahns. In this work, we use a LevenbergMarquardt algorithm to provide this damping. Convergence is further accelerated by a special procedure for generating an initial guess of model parameters. For a GaussNewton or Levenberg-Marquardt algorithm to be efficient, it must converge rapidly and one must employ an efficient algorithm for computing the derivatives (sensitivity coefficients) needed to form the Hessian or modified Hessian matrix. As will be seen, our form of the Levenberg-Marquardt algorithm typically requires on the order of five iterations to obtain convergence and our method for generating sensitivity coefficients is computationally efficient. Landa and Landa and Horne realized that one can directly calculate the sensitivity of pressure to channel parameters, i.e., their procedure avoids calculation of the sensitivity of pressure to gridblock permeabilities and porosities. The appropriate formula relies on an application of the chain rule that allows one to generate the sensitivity of all gridblock pressures to channel parameters by the gradient simulator method (direct method). In their procedure, at each simulator time step, one matrix problem is solved for each channel parameter. The coefficient matrix is the same for each problem, only the right-hand side is different. Thus, if the channel is parameterized using n parameters, we solve a matrix problem with n right-hand sides. If the number of parameters is relatively small, say on the order of a dozen, this should be an efficient procedure. However, the number of geometric parameters used to define our stochastic channel is always greater than four times the number of gridblocks in the x-direction and so it is desirable to find a more efficient procedure for generating sensitivity coefficients. In our work, the procedure of He et al. is applied to calculate the sensitivity of wellbore pressure to gridblock permeabilities and porosities and then these results are inserted into analytical formulas to compute the sensitivity of pressure to model parameters. If Nw is the number of wells at which observed conditioning data has been recorded, this requires solving the reservoir simulator matrix problem with Nw + 1 right-hand sides and then evaluating up to 2N convolution integrals where N is the number of simulator gridblocks. If all conditioning data are measured at the same well, then only one reservoir simulation run is required. This procedure for generating sensitivity coefficients to model parameters and the rapidly converging Levenberg-Marquardt scheme implemented yield an overall procedure for generating realizations that is relatively computationally efficient. Additional information on our procedure for generating sensitivity coefficients is given in Appendix A. In a sense, our work extends some key earlier contributions to the more general problem of conditioning the stochastic channel model of Georgsen and Omre to pressure data using the Bayesian Monte-Carlo approach the last two authors of this paper have championed in recent years. The channel model is more general than the one used by Landa and Horne who considered only twodimensional channels and the one considered by Rahon et al. who assumed the centerline of the channel was a known straight line. Unlike the last two papers, we do not assume the permeability and porosity of the channel and non-channel facies are known a priori. For problems of practical size, the number of model parameters may be on the order of a few hundred, and the gradient simulator SPE 56682 ZHUOXIN BI, D.S. OLIVER AND A.C. REYNOLDS 3 method used by Landa and Horne is not the most efficient choice for generating sensitivity coefficients. Thus, we combine the method of He et al. with analytical expressions to generate sensitivity coefficients; as in the procedure applied by Landa and Horne, the chain rule plays a key role in the development of the formula for sensitivity coefficients. Our set of model parameters includes the geometric parameters which govern the location, size and shape of the channel and channel and non-channel permeability and porosity. Unlike previous work, our formulation is probabilistic; the a posteriori probability density function (pdf) is the conditional pdf for the model given the observed pressure data and well observations of the channel thickness and depth. Moreover, our objective is simulation rather than estimation. By generating multiple realizations of the model, one can evaluate the reduction in uncertainty obtained by conditioning the model to pressure data and/or other information. Realizations of model parameters are obtained with the randomized maximum likelihood method. Each realization requires the minimization of a different objective function which is accomplished by using a special form of the Levenberg-Marquardt algorithm. Stochastic Variables and Probability Density